// Diffractor.h: interface for the Diffractor class.
//
//////////////////////////////////////////////////////////////////////

#if !defined(AFX_DIFFRACTOR_H__B62F0648_4A9C_11D1_93B2_006097737B41__INCLUDED_)
#define AFX_DIFFRACTOR_H__B62F0648_4A9C_11D1_93B2_006097737B41__INCLUDED_

#if _MSC_VER >= 1000
#pragma once
#endif // _MSC_VER >= 1000

#include "DataHolder.h"
#include "Position.h"
#include <cmath>
#include <iostream>
using namespace std;
#include <fstream>
#include <cstdlib>

#include "data_type_to_read_hcv.h"

class Diffractor  
{
public:
	Diffractor(char *DiffractorFile, 
			   Position *W1, Position *W2, Position *W3, char **atomLabel, int NAtom);
	virtual ~Diffractor(void);

	double *Build(Position *AtomPosition);
	int     GetNIntensity(void); 
	void	GetLikelihoodInput(DataHolder* &OutData);

  /// Calculates real and imaginary structure factor values 
  /// and store these in protected att ReAmplitude and ImAmplitude, 
  /// and also populate protected att FDeriv
	virtual void FourierTransform(Position *XValues) = 0;
  /// this method requires FourierTransform has been called first (using FDeriv)
  /// takes as input deriv of chi-sq with respect to intensities and returns
  /// derivs of chi-sq with to orthonormal atom coord
	virtual void TransformDerivatives(double *ID, Position *XD) = 0;


private:
	DataHolder	*Data;
	int			MaxH, MaxK, MaxL; 
	Position	FracX, FracY, FracZ;
	PMetric		M;

	void Diffractor_read_old_style(char *DiffractorFile, 
			   Position *W1, Position *W2, Position *W3, char **atomLabel);
	void Diffractor_read_hcv(char *DiffractorFile, 
			   Position *W1, Position *W2, Position *W3, char **atomLabel);

  // used for reading in f'ing .hcv file
  void get_information_infile(long *number_of_lines_ptr, long *num_extra_lines_ptr, 
							long *highest_multiplet_ptr, char infilename[]);

  FILE* Fopen(char *fname, char *fmode);
  store_data data;
  void populate_store_data(char* DiffractorFile, vector<int>& num_cor_element);
  int count_number_in_line(const char* line);
  
  // transfer form factor decay to data - see description in cpp
  void convertFormFactorDecay(char **atomLabel);

protected:
	int			NIntensity, NAtom, *H,  *K,  *L;
	double	*Intensity;
  // real and imaginary s.f.
  double  *ReAmplitude, *ImAmplitude; 
  double **FiniteFactor;
	double		*SinX, *CosX, *SinY, *CosY, *SinZ, *CosZ;
  // DH, DK, DL are double precision verions of the miller coordinates H, K, L
  // added for speed avoid casting integer to double
	double		*DH, *DK, *DL;
  /// NAtom * 3 matrix. Holds derivatives of chi-sq with respect to fractional coord
  /// As calculated in pure virtual method TransformDerivatives
  double **FDeriv;
	void		Ladder(double *Sin, double *Cos, double S, double C, int N);
	void		SetUpTrig(Position &APosition);
  /// Transform derivatives with respect to crystallographic coordinates,
  /// stored in local variable FDeriv, which is created by TransformDerivatives(), 
  /// to derivatives with respect to orthornormal coordinates times 4*PI
	void		MakeXDerivative(Position *XDerivative);
};


//////////////////////////////////////////////////////////////////////
// inline functions
//////////////////////////////////////////////////////////////////////

inline int Diffractor::GetNIntensity(void) 
{	return (NIntensity);	}


inline void	Diffractor::GetLikelihoodInput(DataHolder* &OutData)
{	OutData = Data;			}


inline void Diffractor::Ladder(double *Sin, double *Cos, double S, double C, int N)
{	
  for(int n=1;n<=N;n++)
  {
		Sin[ n] =  Sin[n-1]*C + Cos[n-1]*S;
		Cos[ n] =  Cos[n-1]*C - Sin[n-1]*S;
		Sin[-n] = -Sin[n];
		Cos[-n] =  Cos[n];
  }
}


inline void Diffractor::SetUpTrig(Position &APosition)
{
double X  = FracX*APosition;
double Y  = FracY*APosition;
double Z  = FracZ*APosition;

Ladder(SinX, CosX, sin(X), cos(X), MaxH);
Ladder(SinY, CosY, sin(Y), cos(Y), MaxK);	
Ladder(SinZ, CosZ, sin(Z), cos(Z), MaxL);
}


inline void Diffractor::MakeXDerivative(Position *XDerivative)
{
for(int n=1;n<=NAtom;n++)
	XDerivative[n] = 2.0*(M*FDeriv[n]);
}



#endif // !defined(AFX_DIFFRACTOR_H__B62F0648_4A9C_11D1_93B2_006097737B41__INCLUDED_)
